%     VIRTUAL HC - программа для получения эффективных свойств сотовых заполнителей
%     Copyright (C) 2018-2020  Михаил Козлов
%     Сайт проекта: <http://www.fastenerdoctor.ru/>
%     Адрес электронной почты для связи с разработчиком:
%     <doctor.fastener@gmail.com>
% 
%     Эта программа является свободным программным обеспечением: Вы можете
%     распространять ее и/или изменять, соблюдая условия Универсальной общественной
%     лицензии GNU, опубликованной Фондом свободного программного обеспечения,
%     либо редакции 3, либо (на Ваше усмотрение) любой редакции, выпущенной позже.
% 
%     Эта программа распространяется в надежде на то, что она окажется полезной,
%     но БЕЗ КАКИХ-ЛИБО ГАРАНТИЙ, включая даже подразумеваемую гарантию ТОВАРНОЙ ПРИГОДНОСТИ,
%     либо ПРИГОДНОСТИ ДЛЯ ОПРЕДЕЛЕННЫХ ЦЕЛЕЙ. Ознакомьтесь с Универсальной общественной
%     лицензией GNU для получения более подробной информации.
% 
%     Вы должны были получить копию Универсальной общественной лицензии GNU
%     вместе с этой программой. Если Вы ее не получили, то перейдите по адресу:
%     <http://www.gnu.org/licenses/>.
% 
%     This program is free software: you can redistribute it and/or modify
%     it under the terms of the GNU General Public License as published by
%     the Free Software Foundation, either version 3 of the License, or
%     (at your option) any later version.
% 
%     This program is distributed in the hope that it will be useful,
%     but WITHOUT ANY WARRANTY; without even the implied warranty of
%     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
%     GNU General Public License for more details.
% 
%     You should have received a copy of the GNU General Public License
%     along with this program.  If not, see <https://www.gnu.org/licenses/>.

clear all; %#ok<CLALL>
clc;
more off;

% #######################################################################################################
% ##                                                                                                   ##
% ##                                        НАСТРОЙКИ ПРОГРАММЫ                                        ##
% ##                                                                                                   ##
% #######################################################################################################

% Для включения возможности вывода информации на русском языке при выполнении программы в среде GNU Octave
% необходимо добавить в файл C:\Octave\Octave-4.4.0\share\octave\site\m\startup\octaverc  следующую строку:
% system("chcp 1251",1);

% ------------------------------ Параметры геометрии ячейки периодичности -------------------------------

% Тип ячейки периодичности
RVE_Type                   = 'TYPE_B'; % TYPE_A | TYPE_B    % Тип ячейки периодичности (см. руководство пользователя)

% Способ задания геометрических параметров ячейки периодичности
RVE_Geom_Model             = 'EDGES_ANGLE'; % EDGE | WRENCH | EDGES_ANGLE | RVE_SIZE | RVE_SIZE_EDGE

        % EDGE - задается размер стороны правильного шестиугольника
        % WRENCH - задается размер "под ключ" правильного шестиугольника
        % EDGES_ANGLE - задаются размеры двух типов сторон шестиугольника и угол раскрытия
        % RVE_SIZE - задаются два габаритных размера прямоугольной ячейки периодичности; длины сторон шестиугольника при этом одинаковые
        % RVE_SIZE_EDGE - задаются два габаритных размера прямоугольной ячейки периодичности и длина стороны шестиугольника с двойной толщиной стенки

        % Если RVE_Geom_Model == EDGE
        Model_EDGE_a = 6; % Размер стороны правильного шестиугольника
        
        % Если RVE_Geom_Model == WRENCH
        Model_WRENCH_h = 3.96875; % Размер "под ключ" правильного шестиугольника
        
        % Если RVE_Geom_Model == EDGES_ANGLE
        Model_EDGES_ANGLE_a = 2.5; % Размер стороны шестиугольника с двойной толщиной стенки
        Model_EDGES_ANGLE_b = 2.5; % Размер стороны шестиугольника с одинарной толщиной стенки
        Model_EDGES_ANGLE_alpha = 120; % [градусы] Угол между двумя соседними сторонами с одинарной толщиной стенки
        
        % Если RVE_Geom_Model == RVE_SIZE
        Model_RVE_SIZE_L = 20; % Размер прямоугольной ячейки периодичности вдоль листов материала (в направлении фольги)
        Model_RVE_SIZE_W = 10; % Размер прямоугольной ячейки периодичности поперек листов материала
        
        % Если RVE_Geom_Model == RVE_SIZE_EDGE
        Model_RVE_SIZE_EDGE_L = 9.33; % Размер прямоугольной ячейки периодичности вдоль листов материала (в направлении фольги)
        Model_RVE_SIZE_EDGE_W = 2.5; % Размер прямоугольной ячейки периодичности поперек листов материала
        Model_RVE_SIZE_EDGE_a = 6.74; % Размер стороны шестиугольника с двойной толщиной стенки

% Способ задания толщины стенки сот
Thickness_Model            = 'FOIL_THICKNESS'; % FOIL_THICKNESS | FROM_DENSITY

        % FOIL_THICKNESS - указывается толщина одинарной стенки сот
        % FROM_DENSITY - толщины вычисляются из объемной плотности сотоблока и плотности материала стенок

        % Если Thickness_Model == FOIL_THICKNESS
        Foil_Thickness = 0.03; % Толщина одинарной стенки сот
        
        % Если Thickness_Model == FROM_DENSITY
        HC_density = 50; % [кг/м3] Объемная плотность сотоблока
        Foil_density = 2700; % [кг/м3] Плотность материала стенок

% Высота ячейки периодичности
RVE_High_Model             = 'AUTO'; % AUTO | SPECIFY

        % AUTO - высота ячейки периодичности вычисляется автоматически
        % SPECIFY - высота ячейки периодичности явно задается пользователем
        
        % Если RVE_High_Model == SPECIFY
        RVE_High_Value     = 15.875;    % Высота ячейки периодичности

% Моделирование наполнителя в ячейках сот
Cell_Filler                = 'FALSE'; % TRUE | FALSE

% ------------------------------------------- Параметры сетки -------------------------------------------

Cell_Mesh_Type             = 'SOLID'; % SHELL | SOLID          % Тип конечных элементов в стенках сот
Element_Order              = 'LINEAR'; % LINEAR | QUADRATIC    % Порядок элементов
N_Elements_Edge            = 15; % Рекомендуемое пользователем число элементов на стороне ячейки сот с одинарной толщиной стенки
N_Elements_High            = 20; % Рекомендуемое пользователем число элементов по высоте ячейки периодичности (используется, если RVE_High_Model == SPECIFY)
Element_Critical_Ratio     = 10; % Максимально допустимое соотношение сторон конечного элемента

% -------------------------------- Параметры массива ячеек периодичности --------------------------------

RVE_Array_L                = 1; % Число ячеек по направлению L
RVE_Array_W                = 1; % Число ячеек по направлению W
RVE_Array_T                = 1; % Число ячеек по направлению T

% ----------------------------------------- Параметры расчета -------------------------------------------

Analysis_Type              = 'MECHANICAL'; % MECHANICAL | HEAT_TRANSFER | MESH_GENERATOR

                           % MECHANICAL - получение эффективных механических свойств с помощью решателя Abaqus
                           % HEAT_TRANSFER - получение эффективных свойств теплопереноса с помощью решателя Abaqus
                           % MESH_GENERATOR - только генерация сетки (без расчета)
                           
Mechanical_BC_TYPE         = 'SOFT'; % SOFT | HARD    % Тип граничных условий в механическом расчете

                           % SOFT - периодические граничные условия со свободной границей ячейки
                           % HARD - периодические граничные условия с фиксированной границей ячейки (использовать с большой осторожностью!)
                           
Constraint_Directions      = 'Z';    % Направления, задающие грани с фиксированной границей в случае Mechanical_BC_TYPE == HARD. Строка, состоящая из символов X, Y, Z.

Chemical_Shrinkage         = 'FALSE'; % TRUE | FALSE    % Расчет эффективных коэффициентов химической усадки помимо основных механических свойств

% ----------------------------------- Свойства материала стенок сот -------------------------------------

Cell_Material_Type         = 'ISOTROPIC'; % ISOTROPIC | ORTHOTROPIC
                           % ISOTROPIC - изотропный материал стенок сот
                           % ORTHOTROPIC - ортотропный материал стенок сот

% ДЛЯ МЕХАНИЧЕСКОГО И ТЕПЛОВОГО РАСЧЕТА

Cell_rho                   = 2.7e-006;    % Плотность

% ДЛЯ МЕХАНИЧЕСКОГО РАСЧЕТА

% Свойства, необходимые в случае изотропных стенок сот
Cell_E                     = 7200;    % Модуль упругости
Cell_nu                    = 0.31;    % Коэффициент Пуассона
Cell_alpha                 = 22e-006;    % Коэффициент линейного теплового расширения
Cell_beta                  = 0;    % Коэффициент химической усадки

% Свойства, необходимые в случае ортотропных стенок сот
Cell_E_1                   = 6600;    % Модуль упругости в направлении 1
Cell_E_2                   = 6600;    % Модуль упругости в направлении 2
Cell_G_12                  = 450;    % Модуль сдвига в плоскости 12
Cell_G_13                  = 400;    % Модуль сдвига в плоскости 13
Cell_G_23                  = 400;    % Модуль сдвига в плоскости 23
Cell_Nu_12                 = 0.03;    % Коэффициент Пуассона 12 (Nu_12*E_2 = Nu_21*E_1)
Cell_alpha_11              = 2.5e-006;    % Коэффициент линейного теплового расширения в направлении 1
Cell_alpha_22              = 2.5e-006;    % Коэффициент линейного теплового расширения в направлении 2
Cell_beta_11               = 0;    % Коэффициент химической усадки в направлении 1
Cell_beta_22               = 0;    % Коэффициент химической усадки в направлении 2

% Дополнительные свойства, необходимые в случае ортотропных стенок сот, моделируемых твердотельными элементами (Cell_Mesh_Type = SOLID)
Cell_E_3                   = 1200;    % Модуль упругости в направлении 3
Cell_Nu_13                 = 0.4;    % Коэффициент Пуассона 13 (Nu_13*E_3 = Nu_31*E_1)
Cell_Nu_23                 = 0.4;    % Коэффициент Пуассона 23 (Nu_23*E_3 = Nu_32*E_2)
Cell_alpha_33              = 70e-006;    % Коэффициент линейного теплового расширения в направлении 3
Cell_beta_33               = 0;    % Коэффициент химической усадки в направлении 3

% ДЛЯ ТЕПЛОВОГО РАСЧЕТА

% Свойства, необходимые для изотропных и ортотропных стенок сот
Cell_Cp                    = 900;    % Удельная теплоемкость

% Свойства, необходимые в случае изотропных стенок сот
Cell_lambda                = 0.235;    % Коэффициент теплопроводности

% Свойства, необходимые в случае ортотропных стенок сот
Cell_lambda_11             = 0.003;    % Коэффициент теплопроводности в направлении 1
Cell_lambda_22             = 0.003;    % Коэффициент теплопроводности в направлении 2
Cell_lambda_33             = 0.0005;    % Коэффициент теплопроводности в направлении 3

% ---------------------------------------- Свойства наполнителя -----------------------------------------

% ДЛЯ МЕХАНИЧЕСКОГО И ТЕПЛОВОГО РАСЧЕТА

Filler_rho                 = 1.2e-006;    % Плотность

% ДЛЯ МЕХАНИЧЕСКОГО РАСЧЕТА

Filler_E                   = 330;    % Модуль упругости
Filler_nu                  = 0.37;    % Коэффициент Пуассона
Filler_alpha               = 35e-006;    % Коэффициент линейного теплового расширения
Filler_beta                = -0.03;    % Коэффициент химической усадки

% ДЛЯ ТЕПЛОВОГО РАСЧЕТА

Filler_Cp                  = 1200;    % Удельная теплоемкость
Filler_lambda              = 0.0002;    % Коэффициент теплопроводности

% % Свойства воздуха для теплового расчета (температура RT, давление 1 бар)
% Filler_rho                 = 1.177E-09;    % Плотность [кг/мм3]
% Filler_Cp                  = 1004.9;    % Удельная теплоемкость [Дж/кг/К]
% Filler_lambda              = 2.624E-05;    % Коэффициент теплопроводности [Вт/мм/К]

% ----------------------------------------- Параметры запуска -------------------------------------------

Work_Directory             = 'M:\Abaqus_Work';    % Рабочая папка для проведения расчета в Abaqus
CPU                        = 4;                   % Число используемых ядер для расчета в Abaqus

% -------------------------------------------------------------------------------------------------------

% #######################################################################################################
% ##                                                                                                   ##
% ##                                     КОНЕЦ НАСТРОЕК ПРОГРАММЫ                                      ##
% ##                                                                                                   ##
% #######################################################################################################

global ErrCode;
global Merge_Tolerance_Multiplier;
global software_version;
global software_date;
global software_year;

ErrCode = 0;
Merge_Tolerance_Multiplier = [];
software_version = [];
software_date = [];
software_year = [];

% Отсчет времени выполнения программы
tic;

% Инициализация базовых переменных
Init_Variables();

% Печать приветствия
Print_Greeting();

% Проверка параметров
switch RVE_Type
    case {'TYPE_A', 'TYPE_B'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, 'Недопустимое значение параметра RVE_Type !\n');
        fprintf (1, 'Разрешенные значения: TYPE_A, TYPE_B\n');
        fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
        toc;
        return;
end

switch RVE_Geom_Model
    case {'EDGE', 'WRENCH', 'EDGES_ANGLE', 'RVE_SIZE', 'RVE_SIZE_EDGE'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, 'Недопустимое значение параметра RVE_Geom_Model !\n');
        fprintf (1, 'Разрешенные значения: EDGE, WRENCH, EDGES_ANGLE, RVE_SIZE, RVE_SIZE_EDGE\n');
        fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
        toc;
        return;
end

if (Model_EDGE_a <= 0)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра Model_EDGE_a !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (Model_WRENCH_h <= 0)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра Model_WRENCH_h !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (Model_EDGES_ANGLE_a <= 0)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра Model_EDGES_ANGLE_a !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (Model_EDGES_ANGLE_b <= 0)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра Model_EDGES_ANGLE_b !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if ((Model_EDGES_ANGLE_alpha <= 0) || (Model_EDGES_ANGLE_alpha >= 360))
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра Model_EDGES_ANGLE_alpha !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (Model_RVE_SIZE_L <= 0)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра Model_RVE_SIZE_L !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (Model_RVE_SIZE_W <= 0)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра Model_RVE_SIZE_W !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (Model_RVE_SIZE_EDGE_L <= 0)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра Model_RVE_SIZE_EDGE_L !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (Model_RVE_SIZE_EDGE_W <= 0)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра Model_RVE_SIZE_EDGE_W !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (Model_RVE_SIZE_EDGE_a <= 0)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра Model_RVE_SIZE_EDGE_a !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

switch Thickness_Model
    case {'FOIL_THICKNESS', 'FROM_DENSITY'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, 'Недопустимое значение параметра Thickness_Model !\n');
        fprintf (1, 'Разрешенные значения: FOIL_THICKNESS, FROM_DENSITY\n');
        fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
        toc;
        return;
end

if (Foil_Thickness <= 0)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра Foil_Thickness !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (HC_density <= 0)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра HC_density !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (Foil_density <= 0)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра Foil_density !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (HC_density / Foil_density >= 1)
    ErrCode = -1;
    fprintf (1, 'Недопустимое соотношение параметров HC_density и Foil_density !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

switch RVE_High_Model
    case {'AUTO', 'SPECIFY'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, 'Недопустимое значение параметра RVE_High_Model !\n');
        fprintf (1, 'Разрешенные значения: AUTO, SPECIFY\n');
        fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
        toc;
        return;
end

if (RVE_High_Value <= 0)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра RVE_High_Value !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

switch Cell_Filler
    case {'TRUE', 'FALSE'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, 'Недопустимое значение параметра Cell_Filler !\n');
        fprintf (1, 'Разрешенные значения: TRUE, FALSE\n');
        fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
        toc;
        return;
end

switch Cell_Mesh_Type
    case {'SHELL', 'SOLID'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, 'Недопустимое значение параметра Cell_Mesh_Type !\n');
        fprintf (1, 'Разрешенные значения: SHELL, SOLID\n');
        fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
        toc;
        return;
end

switch Element_Order
    case {'LINEAR', 'QUADRATIC'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, 'Недопустимое значение параметра Element_Order !\n');
        fprintf (1, 'Разрешенные значения: LINEAR, QUADRATIC\n');
        fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
        toc;
        return;
end

if (N_Elements_Edge ~= floor(N_Elements_Edge))
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра N_Elements_Edge !\n');
    fprintf (1, 'Параметр должен быть целым числом.\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (N_Elements_Edge < 1)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра N_Elements_Edge !\n');
    fprintf (1, 'Параметр должен быть больше или равен единице.\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (N_Elements_High ~= floor(N_Elements_High))
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра N_Elements_High !\n');
    fprintf (1, 'Параметр должен быть целым числом.\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (N_Elements_High < 1)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра N_Elements_High !\n');
    fprintf (1, 'Параметр должен быть больше или равен единице.\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (Element_Critical_Ratio <= 1)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра Element_Critical_Ratio !\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (RVE_Array_L ~= floor(RVE_Array_L))
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра RVE_Array_L !\n');
    fprintf (1, 'Параметр должен быть целым числом.\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (RVE_Array_L < 1)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра RVE_Array_L !\n');
    fprintf (1, 'Параметр должен быть больше или равен единице.\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (RVE_Array_W ~= floor(RVE_Array_W))
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра RVE_Array_W !\n');
    fprintf (1, 'Параметр должен быть целым числом.\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (RVE_Array_W < 1)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра RVE_Array_W !\n');
    fprintf (1, 'Параметр должен быть больше или равен единице.\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (RVE_Array_T ~= floor(RVE_Array_T))
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра RVE_Array_T !\n');
    fprintf (1, 'Параметр должен быть целым числом.\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (RVE_Array_T < 1)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра RVE_Array_T !\n');
    fprintf (1, 'Параметр должен быть больше или равен единице.\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

switch Analysis_Type
    case {'MECHANICAL', 'HEAT_TRANSFER', 'MESH_GENERATOR'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, 'Недопустимое значение параметра Analysis_Type !\n');
        fprintf (1, 'Разрешенные значения: MECHANICAL, HEAT_TRANSFER, MESH_GENERATOR\n');
        fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
        toc;
        return;
end

switch Mechanical_BC_TYPE
    case {'SOFT', 'HARD'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, 'Недопустимое значение параметра Mechanical_BC_TYPE !\n');
        fprintf (1, 'Разрешенные значения: SOFT, HARD\n');
        fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
        toc;
        return;
end

if ((numel(Constraint_Directions) < 1) || (numel(Constraint_Directions) > 3))
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра Constraint_Directions !\n');
    fprintf (1, 'Разрешенные значения: строка, состоящая из символов X, Y, Z.\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

for i=1:numel(Constraint_Directions)
    switch Constraint_Directions(i)
        case {'X', 'Y', 'Z'}
            % Допустимые значения
        otherwise
            ErrCode = -1;
            fprintf (1, 'Недопустимое значение параметра Constraint_Directions !\n');
            fprintf (1, 'Разрешенные значения: строка, состоящая из символов X, Y, Z.\n');
            fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
            toc;
            return;
    end
end

switch Chemical_Shrinkage
    case {'TRUE', 'FALSE'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, 'Недопустимое значение параметра Chemical_Shrinkage !\n');
        fprintf (1, 'Разрешенные значения: TRUE, FALSE\n');
        fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
        toc;
        return;
end

switch Cell_Material_Type
    case {'ISOTROPIC', 'ORTHOTROPIC'}
        % Допустимые значения
    otherwise
        ErrCode = -1;
        fprintf (1, 'Недопустимое значение параметра Cell_Material_Type !\n');
        fprintf (1, 'Разрешенные значения: ISOTROPIC, ORTHOTROPIC\n');
        fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
        toc;
        return;
end

if (CPU ~= floor(CPU))
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра CPU !\n');
    fprintf (1, 'Параметр должен быть целым числом.\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

if (CPU < 1)
    ErrCode = -1;
    fprintf (1, 'Недопустимое значение параметра CPU !\n');
    fprintf (1, 'Параметр должен быть больше или равен единице.\n');
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

% Счетчики числа сообщений
N_Info = 0;
N_Warnings = 0;

% Автоматическое изменение параметров в зависимости от типа расчета
if (strcmp(Analysis_Type, 'HEAT_TRANSFER'))
    
    if (strcmp(Cell_Mesh_Type, 'SHELL'))
        Cell_Mesh_Type = 'SOLID';
        fprintf (1, 'ИНФОРМАЦИЯ: параметр Cell_Mesh_Type был автоматически изменен на SOLID,\n');
        fprintf (1, '            так как будет выполняться расчет свойств теплопереноса.\n\n');
        N_Info = N_Info + 1;
    end
    
end

if (strcmp(Analysis_Type, 'MECHANICAL'))
    
    if (strcmp(RVE_Type, 'TYPE_A') && strcmp(Cell_Filler, 'FALSE'))
        RVE_Type = 'TYPE_B';
        fprintf (1, 'ИНФОРМАЦИЯ: параметр RVE_Type был автоматически изменен на TYPE_B, так как будет\n');
        fprintf (1, '            выполняться расчет механических свойств ячейки без наполнителя.\n\n');
        N_Info = N_Info + 1;
    end
    
    if (strcmp(Cell_Mesh_Type, 'SHELL') && strcmp(Cell_Filler, 'FALSE') && (RVE_Array_W == 1))
        RVE_Array_W = 2;
        fprintf (1, 'ИНФОРМАЦИЯ: параметр RVE_Array_W был автоматически изменен на 2, так как будет\n');
        fprintf (1, '            выполняться расчет механических свойств ячейки без наполнителя,\n');
        fprintf (1, '            в которой стенки сот смоделированы оболочечными элементами.\n\n');
        N_Info = N_Info + 1;
    end
    
    if (strcmp(Cell_Mesh_Type, 'SOLID') && strcmp(Cell_Material_Type, 'ORTHOTROPIC'))
        fprintf (1, 'ПРЕДУПРЕЖДЕНИЕ! Так как параметр Cell_Mesh_Type установлен в значение SOLID, стенки сот\n');
        fprintf (1, '                будут смоделированы твердотельными элементами. Так как для стенок сот\n');
        fprintf (1, '                задан ортотропный материал, то из-за особенностей сетки могут быть получены\n');
        fprintf (1, '                некорректные значения коэффициентов теплового расширения и химической усадки\n');
        fprintf (1, '                по направлениям 1 и 2. В данном случае рекомендуется использовать оболочечные\n');
        fprintf (1, '                элементы для моделирования стенок сот. Для этого параметру Cell_Mesh_Type\n');
        fprintf (1, '                необходимо задать значение SHELL.\n\n');
        N_Warnings = N_Warnings + 1;
    end
    
    if (strcmp(Mechanical_BC_TYPE, 'HARD'))
        fprintf (1, 'ПРЕДУПРЕЖДЕНИЕ! В расчете используются периодические граничные условия с фиксированной границей\n');
        fprintf (1, '                ячейки (Mechanical_BC_TYPE = HARD). Требуется особое внимание при анализе\n');
        fprintf (1, '                полученных эффективных свойств.\n\n');
        N_Warnings = N_Warnings + 1;
    end
    
    if (strcmp(Mechanical_BC_TYPE, 'HARD') && strcmp(Cell_Mesh_Type, 'SOLID') && strcmp(Cell_Filler, 'TRUE') && (~isempty(strfind(Constraint_Directions, 'X'))) && (~isempty(strfind(Constraint_Directions, 'Y')))) %#ok<STREMP>
        fprintf (1, 'ПРЕДУПРЕЖДЕНИЕ! В расчете используются периодические граничные условия с фиксированной границей\n');
        fprintf (1, '                ячейки (зафиксированы грани по X и Y). Так как параметр Cell_Mesh_Type установлен\n');
        fprintf (1, '                в значение SOLID, стенки сот будут смоделированы твердотельными элементами.\n');
        fprintf (1, '                При таком сочетании параметров может быть получено некорректное значение\n');
        fprintf (1, '                эффективного модуля сдвига G12 из-за того, что в углах ячейки принудительно\n');
        fprintf (1, '                создается трансверсальный сдвиг стенок сот. Полученное при этом значение модуля\n');
        fprintf (1, '                сдвига G12 может превышать реальное. Величина этого эффекта зависит от упругих\n');
        fprintf (1, '                свойств наполнителя.\n\n');
        N_Warnings = N_Warnings + 1;
    end
    
    if (strcmp(Mechanical_BC_TYPE, 'HARD') && strcmp(Cell_Mesh_Type, 'SOLID') && strcmp(Cell_Filler, 'FALSE') && (~isempty(strfind(Constraint_Directions, 'X'))) && (~isempty(strfind(Constraint_Directions, 'Y'))))  %#ok<STREMP>
        fprintf (1, 'ПРЕДУПРЕЖДЕНИЕ! В расчете механических свойств ячейки без наполнителя используются периодические\n');
        fprintf (1, '                граничные условия с фиксированной границей ячейки (зафиксированы грани по X и Y).\n');
        fprintf (1, '                Так как параметр Cell_Mesh_Type установлен в значение SOLID, стенки сот будут\n');
        fprintf (1, '                смоделированы твердотельными элементами. При таком сочетании параметров будет\n');
        fprintf (1, '                получено некорректное значение эффективного модуля сдвига G12 из-за того, что\n');
        fprintf (1, '                в углах ячейки принудительно создается трансверсальный сдвиг стенок сот. Полученное\n');
        fprintf (1, '                при этом значение модуля сдвига G12 будет значительно превышать реальное.\n\n');
        N_Warnings = N_Warnings + 1;
    end
    
end

% Необходимость наличия наполнителя в расчете свойств теплопереноса
if (strcmp(Analysis_Type, 'HEAT_TRANSFER'))
    if (strcmp(Cell_Filler, 'FALSE'))
        ErrCode = -1;
        fprintf (1, 'Ошибка: параметр Cell_Filler должен иметь значение TRUE в расчете свойств теплопереноса!\n');
        fprintf (1, 'Если в сотах нет наполнителя, то присвойте ему свойства воздуха.\n');
        fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
        toc;
        return;
    end
end

% Создание базовых точек для оболочечной ячейки
[RVE_Frame_Points, RVE_Frame_Size] = Construct_RVE_Shell(RVE_Type, RVE_Geom_Model, ...
                                                         Model_EDGE_a, ...
                                                         Model_WRENCH_h, ...
                                                         Model_EDGES_ANGLE_a, Model_EDGES_ANGLE_b, Model_EDGES_ANGLE_alpha, ...
                                                         Model_RVE_SIZE_L, Model_RVE_SIZE_W, ...
                                                         Model_RVE_SIZE_EDGE_L, Model_RVE_SIZE_EDGE_W, Model_RVE_SIZE_EDGE_a);

if (ErrCode < 0)
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

% Вычисление толщины одинарной стенки сот
if (strcmp(Thickness_Model, 'FROM_DENSITY'))
    Foil_Thickness = Calculate_Foil_Thickness(Foil_density, HC_density, RVE_Frame_Points, RVE_Frame_Size);
end

% Перестроение базовых точек в случае твердотельной ячейки
if (strcmp(Cell_Mesh_Type, 'SOLID'))
    RVE_Frame_Points = Construct_RVE_Solid(RVE_Frame_Points, Foil_Thickness, RVE_Type);
end

if (ErrCode < 0)
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

% Длина и ширина одинарной ячейки периодичности
RVE_Size = zeros(1,3);
RVE_Size(1) = RVE_Frame_Size(1);
RVE_Size(2) = RVE_Frame_Size(2);

% Вычисление высоты одинарной ячейки периодичности
RVE_Size(3) = Calculate_RVE_High(RVE_High_Model, RVE_High_Value, RVE_Frame_Points);

if (ErrCode < 0)
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

% Объемная доля стенок
V_Cell = Cell_Volume_Fraction(RVE_Type, RVE_Frame_Points, RVE_Size, Cell_Mesh_Type, Foil_Thickness, Thickness_Model, Foil_density, HC_density);

if (ErrCode < 0)
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

% Построение сетки
[NODES, ELEMENTS, RVE_Coordinates, Angle, Node_Merge_Tolerance] = Generate_Mesh(RVE_Type, RVE_Frame_Points, RVE_Size, Cell_Mesh_Type, Element_Order, N_Elements_Edge, Element_Critical_Ratio, Foil_Thickness, Cell_Filler, RVE_Array_L, RVE_Array_W, RVE_Array_T, Analysis_Type, RVE_High_Model, N_Elements_High);

if (ErrCode < 0)
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

% Корректировка Work_Directory
if (Work_Directory(length(Work_Directory)) == '\')
    Work_Directory = Work_Directory(1:length(Work_Directory)-1);
end

% Заполнение структуры со свойствами материала
Material_Properties.Cell_E = Cell_E;
Material_Properties.Cell_nu = Cell_nu;
Material_Properties.Cell_alpha = Cell_alpha;
Material_Properties.Cell_beta = Cell_beta;
Material_Properties.Cell_E_1 = Cell_E_1;
Material_Properties.Cell_E_2 = Cell_E_2;
Material_Properties.Cell_E_3 = Cell_E_3;
Material_Properties.Cell_G_12 = Cell_G_12;
Material_Properties.Cell_G_13 = Cell_G_13;
Material_Properties.Cell_G_23 = Cell_G_23;
Material_Properties.Cell_Nu_12 = Cell_Nu_12;
Material_Properties.Cell_Nu_13 = Cell_Nu_13;
Material_Properties.Cell_Nu_23 = Cell_Nu_23;
Material_Properties.Cell_alpha_11 = Cell_alpha_11;
Material_Properties.Cell_alpha_22 = Cell_alpha_22;
Material_Properties.Cell_alpha_33 = Cell_alpha_33;
Material_Properties.Cell_beta_11 = Cell_beta_11;
Material_Properties.Cell_beta_22 = Cell_beta_22;
Material_Properties.Cell_beta_33 = Cell_beta_33;
Material_Properties.Cell_rho = Cell_rho;
Material_Properties.Cell_Cp = Cell_Cp;
Material_Properties.Cell_lambda = Cell_lambda;
Material_Properties.Cell_lambda_11 = Cell_lambda_11;
Material_Properties.Cell_lambda_22 = Cell_lambda_22;
Material_Properties.Cell_lambda_33 = Cell_lambda_33;
Material_Properties.Filler_E = Filler_E;
Material_Properties.Filler_nu = Filler_nu;
Material_Properties.Filler_alpha = Filler_alpha;
Material_Properties.Filler_beta = Filler_beta;
Material_Properties.Filler_rho = Filler_rho;
Material_Properties.Filler_Cp = Filler_Cp;
Material_Properties.Filler_lambda = Filler_lambda;

% Проверка свойств материала
Check_Material_Properties(Material_Properties);

if (ErrCode < 0)
    fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
    toc;
    return;
end

% Данные для построения систем координат
Coordinate_Systems = Define_Coord_Systems(RVE_Size, Angle);

% Ветвление по типу анализа
switch Analysis_Type
    
    case 'MECHANICAL'
        
        Mechanical_Procedure(Work_Directory, CPU, NODES, ELEMENTS, Cell_Mesh_Type, Cell_Filler, Cell_Material_Type, Material_Properties, Coordinate_Systems, Foil_Thickness, V_Cell, RVE_Coordinates, Node_Merge_Tolerance, Chemical_Shrinkage, Mechanical_BC_TYPE, Constraint_Directions);
        
        if (ErrCode < 0)
            fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
            toc;
            return;
        end
        
    case 'HEAT_TRANSFER'
        
        Heat_Transfer_Procedure(Work_Directory, CPU, NODES, ELEMENTS, Cell_Mesh_Type, Cell_Filler, Cell_Material_Type, Material_Properties, Coordinate_Systems, V_Cell, RVE_Coordinates, Node_Merge_Tolerance);
        
        if (ErrCode < 0)
            fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
            toc;
            return;
        end
        
    case 'MESH_GENERATOR'
        
        Write_Input_Mesh(Work_Directory, NODES, ELEMENTS, Cell_Mesh_Type, Cell_Filler, Cell_Material_Type, Material_Properties, Coordinate_Systems, Foil_Thickness);

        if (ErrCode < 0)
            fprintf (1, 'Работа программы завершена с ошибкой!\n\n');
            toc;
            return;
        end
    
end

% Печать разделителя
Print_Delimiter();

% Сообщение об успешном завершении работы программы
fprintf (1, 'Работа программы завершена успешно!\n\n');

% Вывод числа предупреждений
if (N_Warnings > 0)
    fprintf (1, '\tПолучено предупреждений: %i\n', N_Warnings);
end

% Вывод числа информационных сообщений
if (N_Info > 0)
    fprintf (1, '\tПолучено информационных сообщений: %i\n', N_Info);
end

% Печать новой строки
if ((N_Warnings > 0) || (N_Info > 0))
    fprintf (1, '\n');
end

% Вывод времени выполнения программы
toc;
